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1.  INTRODUCTION 

This  report  is  essentially  an  appendix  to  Ref.  1.  We  shall  show  here  how  Kqs.  (44)  and  (49)  of 
Ref.  1  are  derived  from  Eq.  (43)of  Kef.  1.  The  derivation  isa  lengthy,  purely  mathematical  exer¬ 
cise,  which  requires  no  assumptions  or  approximations  not  already  implicit  in  Eq.  (43)  of  Ref  1 . 

The  article  in  Ref.  1  addresses  the  problem  of  theoretically  calculating  the  instantaneous 
power  hackscuttered  from  a  cloud  in  a  so-called  "monostatic  lidar  system.”  Such  a  system  consists 
of  essentially  two  components:  a  pulse-type  laser,  which  fires  a  short,  highly  collimated  pulse  into 
the  cloud  at  time  t  -  0;  and  a  colocated,  conically  baffled  receiver,  which  measures  backscattered 
radiation  at  times  OO.  The  ultimate  quantity  to  be  calculated  is  </((U),  the  power  measured  bv  the 
receiver  at  any  time  t  due  to  photons  that  have  been  scattered  exactly  n  times  by  the  cloud 
particles,  where  n  is  any  positive  integer.  In  Ref.  1  it  was  shown  that,  for  a  sufficiently  small 
receiver,  J  {/)  is  completely  determined  by  a  certain  function  P  (/,0,0)  Isee  Eqs.  (19),  (21 )-( 23 )  of 
Ref  1 1  It  was  further  shown  that  P((</,0, 0)  has  the  following  formal  representation  |see  Kq.  (43) 
and  Eqs.  (30)-(33)  of  Ref.  1|: 


r  *  r  ft  r  In  r  n  r  In 

du  •••  du  ,  dO  dtp  ■■■  JO  I  i Jtp 

o  0  Jo  "-‘Jo  1  Jo  1  Jo  "Jo  " 

n  ,  t  -  1 

X  cxp\pji\ - (z -e  ri)~  ')|  [  J  /W)sinf)  /(  u^z-e)>b^  H-z-e  >cos </>()) 

i  =  l  ; =0 

x  4  t-c~'  N  u  (1  —  (z  •  e  Hz  •  e  )~  '|  )  A\  ^  u  I (x  ■  e  )  — (z •  e  ) (x  •  e  )(z •  e  )~  1 1  j 

\  *—  i  i  n  /  '  •—  (i  inn / 

i=l)  (=0 

J  "  1 

x  M  ^  u  |(ye  !-(*•«  Kye  )<z-e  (I) 

A  i  i  i  *  n  t\  I 

t  =0 

The  mathematical  meanings  of  the  various  quantities  in  this  formula  are  as  follows:  fi  ,  (},  h  and 
are  any  constants  satisfying 

()</;-S/;.  hid,  0<v0<h/2;  (2a) 


f  is  any  function  satisfying 


!)</! t»  '  r.  for  ()</)<« 


<V  is  the  Dirac  delta  function,  defined  by  the  pair  of  equations 


<V( x  -  xn)  -  0.  if  ,v  r  xn  . 


V  t  )  <V(  x  —  t())  f/.t  S  u  ( ,V()) , 


3 


NWCTP  6605 


for  any  function  g  of  x;  /  is  the  "inequality  function”,  defined  by 


1,  if  "inequality”  is  satisfied, 

/("inequality”)  *  <4) 

0,  if  "inequality”  is  not  satisfied; 


x,  y  and  z  are  the  basis  vectors  of  a  Cartesian  coordinate  frame,  hereinafter  referred  to  as  "the  xyz- 
frame;”  and  finally,  e0,  elf ...,  en  are  unit  vectors  given  by 

eQ  =  z.  (5a) 

e(  =  x(  sin0(cos^(  +  yt  sin0(  a\n<pi  +  z(  cos0( .  (i=l,...,n)  (5b) 


where  the  auxiliary  basis  vectors  x  ,  y  ■  and  z(  are  defined  recursively  according  to 


x,  =  x,  y,  =  y  and  z,  =  z, 
*,  =  e,-i 

y,  =  (zxe(I)/|zxe(,,| 
x  =  y  xz 

i  J  i  i 


(i  =2, ... ,  n) 


(6a) 

(6b) 


Because  of  the  recursive  nature  of  Kqs.  (5)  and  (6),  the  xyz-frame  components  of  e(  will  depend  on 
all  the  angles  0l,<P{,...,Ql,<pi.  A  computer-oriented  procedure  for  calculating  x-e;,  y-e(  and  z-e(  from 
these  angles  is  developed  in  the  Appendix  I  see  Eqs.  ( A7)  |. 

Physically,  the  vector 

u  e(  =  u;  =1TS(  +  ,  (t  =  0,l, ... ,  n)  (7) 


may  be  interpreted  as  the  trajectory  of  an  n-scattered  photon  between  its  ith  scattering,  at  point  S(, 
and  its  (i+  l)th  scattering,  at  point  S(  + ,,  with  S()  and  S  + 1  both  coinciding  with  the  origin  O  of  the 
xyx-frame  [see  Fig.  1 1.  Eq.  (5a)  shows  that  the  photon  initially  leaves  the  origin  O  along  the  +  z 
axis,  and  Eqs.  (5b)  and  (6)  show  that  the  direction  of  u(  for  i  >  1  is  measured  by  the  polar  and 
azimuthal  angles  0t  and  <pt  in  a  frame  whose  polar  axis  points  along  u;  _  , 


Eq.  (1)  is  the  starting  point  for  our  analysis  here.  Simply  stated,  our  goal  is  to  "pacify”  Eq. 

(1 ) - i.e.,  to  reduce  it  to  a  form  that  can  be  evaluated  by  standard  numerical  techniques.  The 

most  obvious  obstacle  to  a  numerical  evaluation  of  Eq.  ( I )  is  the  presence  of  the  three  delta 
functions  in  the  integrand;  these  must  be  analytically  integrated  out.  An  important  tool  for 
accomplishing  these  integrations  is  the  "delta  function  change -of- variable  theorem”  (see  Kef  2  for 

a  proof  of  this  theorem):  If  h  is  a  differentiable  function  ofx  whose  only  zeros  arc  at  x, ,  x., . xm, 

and  if  h'(xr)  *0  for  i=  1, ... ,  m,  then 


<V(/i(x))  = 


<Y(x  — x  ) 

V  _ L 


(8) 
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Before  beginning  the  task  of  eliminating  the  delta  functions,  we  shall  simplify  Kq.  ( 1 ) 
slightly  by  rendering  the  integration  variables  dimensionless.  To  do  this,  we  make  the  scaling 
transformation 

ut~*  ut'  =  ujct.  (j  =  0,...,n  - 1)  (9a) 

Then 

du^ctdu-  (i  =  0,...,n~  1) 


/ 1 -l 

/(  ^  u j  (z  •  e^)> ft  j  =  f(  Uj'iz-  e^)> b/ct  j.  (i=  l ,...,n) 


7=0  7=0 

Furthermore,  since  Kq.  (8)  implies  that  6'(ax)  =  |a|_ltf(x),  then 


11-1 

/  n-l 

(  —  c~ 1  ^  u  a 

l  i  i 

•|  (  —  c~  V(  ^  ul  al )  =  l~  1  < 

>1 

1 

1=0 

1=0 

1  =  0 

,  n  -  I  ,  n  -  t  .  n  - 1 

4  u  t ft  )  =  4  cl  ^  u  'ft  )  =  (c()-IfV|  «  'ft  | 

(=0  (=0  ’  1  (=0 

Substituting  the  above  forms  into  Kq  ( 1 ),  and  then  relabeling  the  integration  variables  u  ’  by 
removing  the  prime, 

(f  =  0 . n  —  1  >  (9b) 

we  obtain 

f  •*>  r  it  2n  ft  r  2n 

rfu0-|  c/a  ,t  _ ,  j  f/0,  (/#,  </fy|  dip  t 

X  cxpf/?ft<  1  -(z-e^)-1)!  [  |  f/KMsint^  u  (z  •  e^)>ft/r(  j  /( -z •  e^ >cos<p0) 

i=l  7=0  7 

/  'I".1  / 

Xml  —  ^  u  |l—  (z  •  e  )  (z  •  e  )~ 1 | ) ftl  ^  u  |(x  •  e  )  —  (z  •  e  )(x  •  e  )(z-e  )~'l| 

'  — ■  (  t  n  /  \  —  i  i  i  fi  n  / 

1=0  1=0 

/  "  .' 

X  U(  |(y  •  e  )-(z •  e  )(y  ■  e  Hz  •  e  ) ~  '|  j .  (10) 

1=0 

For  brevity  we  shall  henceforth  refer  to  the  three  delta  functions  in  Kq  (10)  as.  reading  from 
left  to  right,  the  (  delta  function,  the  x  delta  function  and  the  v  delta  function  |cf  Kq.  (43)  of  Kef  1 1 
We  turn  now  to  the  task  of  integrating  out  these  della  functions 
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2.  ELIMINATING  THE  x-  AND  y-DELTA  FUNCTIONS 
We  define  the  vector  variables  A, ,  A2,...,  Art  by 

- 1 

(i=  !,...,«) 


i-l  i-i 

A  =  ^  u  =  N  u  e 


i  j  j  j 

7=0  J= 0 


(11 


Geometrically,  A .  is  the  (scaled)  position  vector  OSi  of  the  point  at  which  the  scattering  occurs 
fsee  Fig.  1 1.  We  note  in  particular  that  A/(  is  independent  of  en,  and  hence  of  the  two  integration 
variables  0n  and  <pn.  In  terms  of  the  vectors  A  ,  Eq.  (10)  can  be  written 


r  «  r  *  r  n  r  In  r  n  r  In 

>  (t,0,0)  =  pnc(ct)n-3exp(-Pct)  \  </«-•-  du  ,  d0t\  dtp  dd  ,  d<p  , 

"  *  Jo  0  Jo  "-'Jo  'Jo  '  Jo  "-'Jo  "-' 

(1-1 

II 


1  =  l 


)  sin()  HA  >b/ct )  L, 

i  i  ij. 


(12) 


where  we  have  defined 


fn  r  'in 

dd  dtp  exp[pbt\  —  (z -e  ) —  1  )|  /f0  )sin0 

0  "Jo  "  "  " 

,  K  - 1 

X  HA  >b/ct)I(-z-e  >cosu>„)  6\  1  -  ^  u  +A  (z-e  )“') 

njt  n  M)  \  i  nj.  n  / 

i=0 

xd\A  —A  (x  •  e  )  (z  •  e  )~ ' )  m  A  -A  (y  •  e  )(z  •  e  )~  ' ).  (13) 

\  njc  n ^  n  n  J  \  n%y  n*z  n  n  f 

The  quantity  L  in  Eq.  (13)  is  seen  to  be  an  integral  over  all  directions  of  the  vector  er).  Our 
goal  in  this  section  is  to  evaluate  L  analytically,  and  in  the  process  eliminate  the  x-and  y-delta 
functions.  To  this  end,  we  first  change  the  integration  variables  in  L  from  (0, ,.$„),  the  polar  and 
azimuthal  angles  of  in  the  x^z^-frame  [see  Kqs.  (5)  and  (6) |,  to  (0,0),  the  polar  and  azimuthal 
angles  of  ej(  in  the  xyz-frame.  This  change  of  variables  allows  the  dot  products  involving  e;(  in  the 
integrand  to  be  written  in  the  relatively  simple  explicit  forms 

x  •  en  =  sinflcosd>,  y  •  e;(  =  sinf)  sin#,  z  •  e  =  cost).  (14a) 

The  Jacobian  of  the  transformation  (0  ,<p  )~*U),tp)  is  such  that  the  form  of  the  differential  solid 

nT  n  r 

angle  is  preserved,  so  we  have 


W\n0 nd0 ndtp n  -  sin OdOdtp. 


(14b) 
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Since  On  appears  as  an  argument  of  f  in  Eq.  (13),  we  will  also  need  a  formula  for  0n  in  terms  of  0 
and  <p.  Such  a  formula  can  be  obtained  by  first  noting  from  Eqs.  (5b)  and  (6)  that  cos 0n  =  z  •  en  = 
ej(  _  .•  e  ;  expanding  the  latter  dot  product  in  terms  of  xyz-components,  we  obtain 

0n  =  O'ld,  </>)  =  arccosl(x-efl  j)sin0cosd>  +  (y -e^jlsinOsin#  +  (z  •  en  _  t )  cos0|.  (15) 

Now  substituting  Eqs.  (14)  into  Eq.  (13),  we  get 


/,  =  |  dO  j  d0exp|/?6(l  — sec0)|/(0  ((0,$))sin0/(A(  ^>fi/c0/(  — cos0>cos(po) 


if  1  —  ^  u  +  A  sec0  I  if  A  — 

\  •fc—  i  nj  )  \  njc 


A  tan0cosd> 

'ix 


A  ^tanO  sin</>  j, 


in  which  0j((0,  (p)  is  understood  to  be  the  function  in  Eq.  (15).  The  dependence  of  the  integrand  on 
the  integration  variables  0  and  <p  has  now  been  rendered  completely  explicit,  so  we  can  perform 
these  two  integrations.  We  shall  integrate  first  over  </>  with  the  help  of  the  y -delta  function,  and 
then  over  0  with  the  help  of  the  x -delta  function. 


Since  -cos0  =  cos(n  -0),  the  second  /-function  in  Eq.  (16)  essentially  requires  that  6>n- 
This  requirement,  co  'pled  with  the  fact  that  <p0<«/2  means  that0>rc/2,  so  we  can  increase  the 
lower  limit  on  the  0-integration  from  0  to  n/2.  We  rewrite  Eq.  ( 16)  in  the  iterated  form 


•it  ,  n  - 1 

L  =  d0  expl/iWl  —  sec0)|sin0  /(A ^  >b!ct)  /(  —  cos0>cosjp())<^  1  —  u(  +  A  ^secO  j //,  (17a) 


where 


'In 

L’  =  d<f>  flO  (0,<h))fHA  —A  tan0 cosd»)tf(A  —A  tan0sind>). 

(I  n  "J  "•*  ',V>'  'IX 


The  details  of  evaluating  first  L‘  and  then  L  turn  out  to  depend  upon  the  signs  of  Anjc  and 
A  „v,  i.e,  upon  which  quadrant  of  the  xy-plane  the  vector  ()S;i  projects  |see  Pig.  l|.  Butitisclear 
from  symmetry  considerations  that  the  final  result  must  be  the  same  in  all  cases.  Therefore,  we 
carry  out  here  the  analysis  only  for  the  case  in  which  Sf|  lies  above  the  first  quadrant  —  i.e. ,  the 

Case:  A  >0  and  A  >0  (18) 

'l.JT  'I.V 

Denoting  the  argument  of  the  y-delta  function  in  Eq  (17b)  by 

hi(p)  -  An  y-  An  .ti*n0sind>,  (19) 

we  rewrite  that  equation  as 
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f2" 

//=  j  d0/W^<0))dM^-A^tan0costfd/(|A(J<|A  ^tanfl|)tf(/i(^t),  (20) 

where  the  /-function  makes  explicit  a  condition  that  is  clearly  necessary  if  h  is  to  vanish  for  some 
value  of  ip.  Assuming  that  condition  is  satisfied,  it  follows  from  Eq.  (19)  that  there  will  in  general 
he  two  angles,  <p{i)  and  $<2>,  for  which  h  vanishes,  namely  those  angles  for  which 

sin<A,l)  =  sin$(2!  =  A  /(A  tan0).  (21a) 

r  r  n.y  n  ,z 

Since  A„^>0and  6t(nl 2,  n )  [cf.  Eq.  (17a)|,  while  A„_v>0  [cf.  Eq.  ( 18)  |,  the  quantity  on  the  right 
side  of  Eq.  (21a)  is  negative.  This  means  that  <p{  l>  and  <£<2'  must  lie  in  the  third  and  fourth 
quadrants,  respectively;  thus,  using  cos 2<p'l>  =  1  -sin2^",  we  deduce  that 

cos<il"  =  (A  2tan20  — A  2  )1/2/(A  tanO)  =  —cos ip'"'  (21h) 

n,z  n  ,v  nj 

Mow,  from  Eq.  (19)  we  have  h'(<p)—  -  A,,^  tan0cos</>,  so  using  Eq.  (21b)  we  deduce  that 

|/T<0"’)|  =  \h'{<pi2,)\  =  (A  (2  tan20-A;i2  )m. 

Therefore,  we  may  use  the  rule  in  Eq.  (8)  to  write 

2 

<5 \h(<p))  =  (A  2  tan20- A  2)'12  V  dUp-ip'").  (22) 

nj  nj  — 

i  =  l 

Substituting  Eq.  (22)  into  Eq.  (20),  the  <p  -  integration  becomes  trivial  (cf.  Eq  <3b)|:  it  yields 

2 

L’  =  /(I A  |<| A  tanfl|)(A  2tan2«-A  2T,r2  V /10  <0,0''’)>A<A  -A  tanflcos*/.1'’).  (23) 

1  ny’  1  n.z  n ^  nK\  *—  n  nji  "j 

i  =  1 

When  the  expressions  for  cos <ptl>  in  Eq.  (21b)  are  substituted  into  the  delta  functions  in  Eq.  (23)  we 
find  that,  for  the  case  A„  f  >0  being  considered  here,  the  argument  of  the  t  =  2  delta  function 
never  vanishes.  That  term  may  therefore  be  dropped,  and  Eq.  (23)  reduces  to 

L'=U\  A  |<|A  tanO|)(A  2  tan20  — A  2  r,/2/10  (O,0"’)><V<A  -|  A  2  tan20  -  A  2  |l/2),  (24) 

1  ft ,  v  n,z  1  njz  nx\  n  n.x  nj  n.x 

where  0t  «),<p' 1  ’)  is  found  from  Eqs  (15)  and  (21 )  to  be 

0  <0,d/ 1  )  =  arccos|(x  -  e  , ) cosd  A1  ( A  “tan“0  — A  )1- 

n  n  -  I  n^z  >i.z  n,\ 

+  (y  •  e  ,)cosWA-,A  +(z-e  ,)cos0|.  (25) 

>1-1  <  -  1 

Substituting  the  above  integrated  form  of  /,'  into  Eq  ( 1 7  a  I ,  we  obtain 
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J^|A  ta  n0|>a(i  -  X  «,  +  *„,  ..sect)  j 
1,2  !=() 

X  expl/tbd  —scc0)l/fO  (6,<ptl)))sin9(A  2  tan 20  — A  ~ )~  1/2  S(hid)) ,  (26) 

'  n  ,T  n,y 

where,  in  anticipation  of  using  Eq.  (8)  again,  we  have  redefined  h  to  be 

M0)  =  A  -(A  2tan20-A  2),/2.  (27) 

Noting  that  the  integration  variable  6  is  confined  to  the  second  quadrant,  where  the  tangent  is 
negative,  we  see  that  h  vanishes  only  at  the  angle  (t 1 1  defined  by 

tant/1’  =  —(A  2  +  A  2  )1/2/A  (28a) 

hjl  n.y  n<z 

This  fact,  coupled  with  the  second  quadrant  restriction,  implies  that 

sint)"’  =  (A  2  +  A  2  )l/2/A  ,  cost/  =  -  A  /A  .  (28b) 

nj  n%y  n  n  n 

The  above  formulas  in  turn  imply  that 

sect)"’  =  -  A  /A  ,  | A  tant/  *'|  =  (A  2  +  A  2  )l/2,  (A  2  tan2(/"-A  2  ),/2  =  A  .  (28c) 

n  >u  >;jr  >\r  u  nj  n.x 

Now,  the  derivative  of  the  function  in  Eq  (27)  is 


•  n 

L  -  HA  >b!ct)  <#)/(- cost)  >  cosip ,.)/(|A 

i  U  „  0  n 

nl 2 


h'(0)  -  —  \ (A  "tan'O  — A  2  )  1/22A  2tan0sec2f) 
Putting  0  =  (/ 1 1  and  using  Eqs.  (28),  we  find  that 


|/i’<f/")|  =  (A  2  +  A  2  )I/2A  2 A  A 

»>  »*,v  ’i  njc  n,z 


The  rule  in  Eq  (8)  then  allows  us  to  write 


ti(hm  -  A  A  A~  2  (A  2  +  A  2  )  -  md(0 -(/ ") . 

njc  njs  n  n.x  n.y 


(29) 


Substituting  Eq  (29)  into  Eq.  (26),  the  0-  integration  is  now  trivially  performed,  with  the  result 


"  -  l 

I.  =  HA  >b/ct)H  —  cost)' 1  '>cos7'  ) /(|A  |<| A  tant/ 1  ’j) 1  —  ^  «  +A  sect/1') 

n.2  T0  >i,v  a  j  \  —  i  n  z  ) 


l  =0 


ex  pi  phi  I  —  sect)  1 ')!/(().  (f)' 1 '  ,<p  1 '))  sint)  1  (A  f  “la  n't)' 1 '  —  A  )  ,r“ 


x-  A  A  A  '-(A  ’  +  A  2)  12 

n  ,»  n  '(.r  *t,v 


(30) 
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Substituting  for  & 1 '  from  Eqs.  (28)  and  simplifying,  we  finally  obtain 

.  "  - 1 

L  =  I(A  > bid)  HA  /A  >cos^Jexp(/J6(l  +A  /A  )|/T0  )A  A-3*!  I  -  v  u  -  A  ),  (31) 

nj  n,z  n  0  r  1  n  n ^  a  n,z  n  \  —  i  >i ) 

1=0 

where  6n  is  now  given  by 

Qn  =  arccos(  -  en_  ,•  AJAn) .  (32) 

Notice  that  the  result  in  Eqs.  (31)  and  (32)  is  manifestly  independent  of  the  signs  of  A,lj:and 
A„tV,  and  hence  of  the  case  assumption  in  Eq.  ( 1 8);  indeed,  if  the  foregoing  analysis  is  repeated  for 
the  cases  in  which  either  or  both  of  AltJC  and  A„  v  are  negative,  the  final  result  will  still  be  Eqs.  (31 ) 
and  (32). 

Substituting  Eq.  (31)  into  Eq.  ( 12),  and  also  putting  sin ftdO  =  -dcosfl  for  i  =  1  -  1,  we 

obtain  the  following  expression  for  Pi  1,0,0): 

r  *  I  -n 

(/«  ••  du  dcosG  d<p  ■■■ 

0  0  j  0  - 1  Jo 

(■  1  rL2rr  .  n  -  1 

X  dcosO  ,  dtp  .  (  I  I  /(A  >b/ct)\liA  >b/ct)HA  /A  >cos<f. .) 

J  J  rl-  1]0  *>l-l\  ■  I  l,Z  )  'I.Z  n  T0 

/  n  >;-l 

x(  [|/(0|))expl/iWH-Afi/A^)|  A^A;3^  1  -  ^  a,  -  A  J.  (33) 

1=1  ”  i=0 

In  this  equation,  the  vectors  A,,..., A  are  given  by  Kq.  (11),  and  0);  is  given  by  Eq.  (32). 


3.  ELIMINATING  THE  t-DELTA  FUNCTION 

The  task  of  integrating  out  the  remaining  delta  function  [  formerly  the  delta  function 
involving  t  in  Eq.  (1)1  is  accomplished  differently  for  n  —  1  than  for  n  22.  We  consider  first  the 
relatively  simple  n=  1  case,  for  which  Eq.  (33)  reads 

r 

P^t,  0,0)  =  fig  c(rt)_2exp(  —  /Set)  du^HA^  T>b/ct)HA^/A{>  cost/.(() 

x  /(O^expl/lW  1  +  A ( /  A  (  )|  A  (  A  “ '* *W  1  —  //()—  A  () 

Since,  according  to  Eqs.  (11)  and  (5a), 
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A  I  Moeo  -  u<)7'  ’ 

then  A\^  =  A\  =  Uq,  and  Kq.  (32)  gives 

0,  =  arccos(  -  eQ  •  uHe0/u())  =  arccos(  -  1)  =  n  . 

The  condition  >cosy>(,  isjust  1  >cosy;0>  which  is  always  satisfied.  Kq.  (34)  therefore 

simplifies  to 

P((/,0,0)  =  /(  cict)~~ cxp(-pet)  |  du0HuQ>b/ct) /frc)exp(2/?6)  a  “  2  5(1  —  2u()) . 


Using  the  rule  in  Kq.  (8),  we  have 


Therefore,  the  above  equation  is 


5(1  -2u0)  =  2-1 5(u0-  1/2) 


P/Z.O.O)  =  2-l/3.c(c<)"2exp|-/i(c/-26)IA«) 


du.Kun>b/ct)  un  25(m„—  1/2  ) 


The  «0-integration  is  now  trivially  performed  with  the  aid  of  Kq.  (3b),  giving 


P  ,(/,<),<))  =  Hct>2b)2cli  (ct)  2exp|  —  ftict— 26)|  fin) .  (35) 

This  result  agrees  exactly  with  Kq.  (28)  of  Ref.  1,  which  was  obtained  from  comparatively  simple 
physical  arguments.  This  agreement  constitutes  a  reassuring,  if  somewhat  limited,  consistency 
check  on  our  calculations  thus  far. 

We  now  turn  to  the  more  interesting  and  challenging  case  in  which  2.  As  a  prelude  to 
eliminating  the  remaining  delta  function  in  Kq.  (33)  for  that  case,  we  make  two  more  integration 
variable  transformations.  The  first  of  these  is  the  transformation 

{un-\<  cos0„  - 1 ’  -  l ]  {A«’  cosV.  D,  (36a) 

where  ip  and  /?  are  the  polar  and  azimuthal  angles  of  the  vector  \ri  in  the  xyz-frame.  Since  both 
sets  of  integration  variables  in  Kq.  (36a)  simply  integrate  the  point  Stover  all  space,  we  have  the 
differential  relation 

u  2  ,du  ,  </cosb  < Id)  =  A  1  dA  dcosi/’t/q.  (36b) 

Using  this  relation  and  the  fact  that  A„  A„  cosy>,  Kq  (33)  becomes,  for  n> 2, 


SSA  giSttJ-J’SJ 
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f'/>  r  '*)  r  00 

du  •  ••  da  dA 

o  0  Jo  "-’Jo  " 

fl  r  2/f  r  1  r  2n  r  l  r  2rc 

_1‘faoH0tJo  J  _lrfCOS0n-2|„  ^„_2|  j^’Jo  dn 

,  n-\ 

x|  j  |  HA  ^>b/ct)^l[Ancosip>b/ct)Hcostf)>cosqi0) 


=(fiH 


expI/Jfe(l  +  secy41  cos</> 


(  "  ~ 1 

s <I’U~*X  <s(  1  -  v  ut  (n> 2)  (37) 


in  which  it  is  understood  that  the  integration  variables  cos 0(  and  <pt  are  absent  if  n  =  2..  The 
quantities  appearing  in  the  integrand  are  given  in  terms  of  the  integration  variables  as  follows. 
The  (n-  1)  vectors  ©0, are  calculated  from  the  recursion  relation  in  Kqs.  (A7).  The(n-l) 
vectors  A,,...,An_,  are  defined  through  Eq  ( 1 1): 

<-l 

A  =  ^  ue  ■  (i  =  l,...,n  —  1)  (38a) 

<  —  j  j 

j  =  n 

The  vector  u^  _  p  and  its  associated  magnitude  un  t  and  unit  direction  vector  e„  _  , ,  are  also 
defined  through  Eq.  (11),  but  now  written  for  i  =  n  in  the  form 

n-2 

u  ,e  ,  =  u  =  A  —  ^  «  e  ,  (38b) 

n  -  1  ft  —  1  rc  -  1  *—  j  i 

J  = ') 

And  finally,  6n_ ,  and  0n  are  given  by  Isee  Eqs.  (5),  (6)  and  (32)1 

0  .  =  arccos(e  •  e  .),  (38c) 

n  -  1  n  -  1  n  -  1  * 

0  =  arccos(-e  .’A  /A  ).  (38d) 

n  n  —  1  n  n 

Our  second  transformation  of  integration  variables  is  another  scaling  transformation, 


which  evidently  gives 


it  v  =  u  IA  ,  <i  =0,...,n -2) 

f  i  in 1 


dn..  -du  .,dA  =  A  "  1  dvn  ■■■  dv  ..dA  . 

0  n  - 1  n  n  0  n  -  2  n 


This  transformation  has  the  effect  of  making  Ati  the  "unit  of  length”  for  all  distance  vectors.  Thus, 
for  i  =  0,...,  ft -2,  we  have  u  -  A  v  ,  where 

v  -ve  it  ~  0,...,  n  -  2)  (40a) 

i  it  * 

And  using  Eq.  (38b)  we  see  that  we  can  also  write  un  j  -  A^vfj  j ,  provided  we  define  v  x  by 


•  ■  •  *  •  *  »  ■  «  *  p  •  ,  *  •  ■  •  *  »'•*#*,  1  »  s  .  *  ,  *  •  *  ^  *  .  *  ,  *  ,*•  t  *  k  .  »  .  • . 

•  •  ♦  •  •  •  - . *  •  *  •  •  •'  •*  -  •  v*  -  *  *  •  *•  •  ■  •  *  *  •  •  •  *  •  *  •  •  *  •  *  -  *  v- • .  * .  *.• 
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n- 2 


v  ,e  ,  =  V  ,  =  a  -  >  e  e  , 

n  - 1  ii  - 1  ii  - 1  — •  j  j 


7  =  0 


where  a  is  the  unit  vector  in  the  direction  of  A  : 


a  =  A  /A  =  xsinipcosq  +  ysimpsinq  -f  z  cosip. 


Finally,  we  see  from  Eq.  (38a)  that  we  can  write  A(  =  AnB;  for  i  =  1,...,  n-  1,  where 


i-l  i-l 

=  V  „  =  V 


B  =  N  v  =  N  ue  .  (i=  1) 

'  —  j  <—  j  j 

1-0  7=» 


With  Eqs.  (39)  -(42),  the  integral  in  Eq.  (37)  takes  the  form 


P  ((,0,0)  =  I) "  c(ct)n  ^exp  (—[Set)  [  dv.-  i  dv  |  A'1  *dA 

i«  ''s  r  H  o  Jo  '1  -  -  ]  0  n 


fl  in  1  rStl  r  1  OX 

dco s0(  drp]  </cosO(i  2J  dip  .,  I  denstf)  dq 

-l  0  J  -  l  •)()  "J  ciisiii,,  0 


'in 


ci>si/i0 


ii  -  1 


X ^  |  |  /(A  B  >h/c()  |  /(A  cosy  >  h/rt) 
i  =  l 

x(  |  |  /10)  jexpl/Ml  -t-secy)!  cosy  A  “  1  u  ~ _“j  «5i(  1  -A  ^ 


i  =  i 


n  -  1 

1  +  V 


1  =  0 


(40h) 


(411 


(42) 


j.  (n> 


2)  (43) 


Now  we  define 

n-  i 

Vsi  +  v  t, 
1  =  0 


(44) 


Then  the  delta  function  in  Eq  (43)  can  he  written 

u  - 1 

/  i 

ft 

\  n\ 

i=0 


>{  1  —  A  1  +  V  v  I  )  =  <K1  -  VA  )  =  V~'ttA  -  V-1), 

\  n  — r  '  I  /  **  *i 


(45) 


where  the  last  step  follows  from  the  rule  in  Eq.  (8)  When  Eq.  (45)  is  substituted  into  Eq.  (43),  the 
Af  integration  can  be  trivially  accomplished:  the  delta  function  is  thereby  eliminated,  and  A  is 
everywhere  replaced  by  V  1 


Before  writing  down  the  result  of  the  A  -integration,  we  want  to  do  two  more  things  to  Eq. 
(43).  The  first  is  simply  to  replace  i/cosy  by  -sinydy  The  second  is  to  fix  the  orientation  of  the  xa 
plane,  which  thus  far  has  not  been  specified.  Owing  to  the  symmetry  of  our  problem  about  the  z- 


uNlV-N  V. ' 


f*  V . 


£££& t. 
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axis,  which  exists  because  the  laser  and  receiver  both  lie  on  that  axis,  we  can  choose  the 
orientation  of  the  xz-plane  freely  without  affecting  the  integrand.  Let  us  now  stipulate  that  the  xz 
plane  contains  the  vector  a.  This  implies,  firstly,  that  ri  =  0  in  Eq.  (41),  and  secondly,  that  the  q- 
integration  in  Eq.  (43)  can  be  replaced  by  a  simple  factor  of  2 n. 

Performing  all  the  operations  described  above,  Eq.  (43)  becomes 

(■*>  i>«  ft  f  2n 

di;0  -|  dai2|  dcosOj  J  d<P{  - 

f  1  f hl  [ /  " _  1  > 

X  dcos0(i  J  d<t>n  dy  (  [  |  nV~lB.z>bfct)jnV~leosip>h/cn 

- 1  0  O  I  =  I 

/  11 

Xexp(/?ft(l +scce/))|  (  |  |  /W  )  jeosysiny’  V'-'"-"1  o~^  ,  (n>2)  (46 

i  =  i 

in  which  it  is  understood  that,  for  n  =  2,  the  integration  variables  cos0(  and  are  absent.  The 
quantities  in  the  integrand  of  Eq.  (46)  arc  related  to  the  integration  variables  according  to  the 
following  specifications:  The  unit  vectors  e0, ...,  en  2  are  found  from  the  recursion  relation  in 
Eqs.  (A7),  while  the  quantities  en  _  J(  a  (,  (9  j,  6n,  Bj, ...,  B/t_ ,  and  V  are  determined  from  the 
formulas: 


a  =  x  sini p  +  z  cosip; 

(47a) 

n-2 

v  ,e  ,  =  a  —  ^  v  e  ; 

n-  1  n  —  l  j  j 

J  =  o 

(47b) 

°„-l  =  arccos(e„-2,e-,-l); 

(47c) 

0  =  arccos(  -  e  ,  •  a); 
n  n  -  l 

(47d) 

£-1 

B  =  ^  v  e  :  U~  1 1 ) 

(47e) 

<  —ii 
i=<) 


•i-i 

V  =  I  +  >'  v  . 


(470 


F  =0 


The  content  of  the  above  relations  is  summarized  geometrically  for  the  cases  n  =  2,  3  and  4  by  the 
diagrams  in  Fig  2 
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4.  BOUNDING  THE  INTEGRAND  ANI)  INTEGRATION  DOMAIN 

In  Eq.  (46)  we  have,  at  last,  an  expression  for  P  it,  0,0)  for  n  22  that  is  free  of  delta  functions. 
It  is  evidently  a  (3  n  -  4)-dimensionul  integral,  the  complexity  of  which  will  usually  dictate  that  it 
be  evaluated  numerically  rather  than  analytically.  However,  Eq.  (46)  is  not  suitable  for  numerical 
evaluation  for  two  reasons:  First,  the  integration  domain  is  unbounded,  since  the  integration 
variables  i>(),...,  vn  _2  have  infinite  ranges;  and  second,  the  integrand  is  unbounded,  since  un_  (  can 
vanish.  ISince  sec<p-*»  as  tp-*n/2,  the  exponential  factor  involving  secip  in  the  integrand  might 
also  seem  to  present  a  boundedness  problem;  however,  the  last  /-function  in  the  integrand  imposes 
the  condition  bsectp<ct/V  <ct,  where  the  last  inequality  follows  from  the  definition  of  V,  so  the 
exponential  is  not  a  problem  !  It  must  be  emphasized  that  these  unbounded  features  of  Eq.  (46)  do 
not  imply  that  the  expression  therein  is  mathematically  ill-defined.  But  they  do  imply  that,  if  we 
want  to  evaluate  the  integral  using  conventional  numerical  techniques,  we  will  have  to  subject  it 
to  some  integration  variable  transformations  that  render  the  integrand  and  the  integration 
domain  bounded. 

The  key  to  obtaining  a  set  of  integration  variables  for  which  the  integrand  and  integration 
domain  are  bounded  turns  out  to  be  the  set  of  vectors  C0,  C,,...,  C(i  ,,  where  C(  is  defined  to  be  the 
vector  from  point  S(  to  point  S([  (see  Fig.  3): 

fa.  if  i  =  0, 


i  - 1 

^  v  e  ,  if  i  =  1  ,...,n  —  1 


Notice  in  particular  that 


C0  =  I  and  C(  _ , 


v  ,  e  . . 

n  -  1  h  -  I 


These  vectors  C(  will  not  themselves  be  our  new  integration  variables,  but  they  will  be  crucial  for 
defining  those  new  variables.  Essentially  what  we  are  going  to  do  now  is,  first,  replace  each  pair  of 
integration  variables  (0(,<£()  in  Eq.  (46)  by  a  new  pair  of  variables  (O' and  second,  replace  each 
length  integration  variable  by  a  new  angular  variable  v( 

The  variables  0  and  <f»i  were  defined  to  be  the  polar  and  azimuthal  angles  respectively  of  the 
vector  e(  relative  to  a  coordinate  frame  whose  polar  axis  is  z(  =  e(  (.  We  now  change  integration 
variables  according  to 

(cost) ,  fb  I  -*<cos()’  ,  (b‘  ),  (i  —  I - 2)  (50a) 

I  r '  l  ~  l 
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where  B't  and  tp' { are  the  polar  and  azimuthal  angles  respectively  of  e(  relative  to  a  coordinate 
frame  whose  polar  axis  is  z'(  =  Ct  1C ..  The  Jacobian  of  this  transformation  is  such  that  each 
differential  solid  angle  element  dcosO  dtp  is  preserved,  so  we  have  the  differential  relation 


(50b) 


dcosSj  d<pl  -dcosOn  2dtpn  2  -  dcos$'}  d#',—  dcos0'n2  ^  „_2 

=  (sin0',— sin0'ri  2)  d0'  x  d<P' {-d0' n_2  dtp'  n  _2 


Since  the  integration  limits  on  0  and  <pt  encompass  all  possible  directions  of  e(,  it  follows  that  O' t 
and  <p'  •  will  have  the  same  respective  limits.  By  definition,  O' ( is  the  angle  between  e(  and  C(  for 
i=  l,...,  n-2,  and  we  can  evidently  extend  that  definition  to  i  =  0  by  simply  defining 


0'()  =  ip  (51) 

The  geometrical  relations  between  the  old  and  new  polar  angles  are  illustrated  in  Fig.  3.  The 
orientation  of  the  azimuthal  plane  that  defines  the  zero  of  tp'  ( is  open,  but  there  will  be  a  minimum 
of  computational  work  later  on  if  we  take  this  plane  to  he  the  one  defined  by  C;  and  z.  Thus,  we  are 
essentially  transforming  from  the  xp^-frame  of  Eqs.  (5)  and  (6)  to  the  x'^'^-fraine  defined  by 


z'  =  C./C, 

y'(  =  (zxC  )/|zxCJ 

*'l  =  y'lxz'r 


(i-  1 . n-2) 


<52a) 


relative  to  which  e(  has  the  component  representation 

e  =  x'  sinfl'  cos#'  4-  y'  sin  O'  sin  tp'  +  z'  cost)'  .  (i  =  1, ... ,  n  -  2)  (52b) 

A  detailed  procedure  for  calculating  the  xy’-components  of  e(  from  the  xvz-components  of  C(  and 
the  angles  0';  and  tp' f  is  developed  in  the  Appendix  |see  Kqs  ( A8)| 


/Mt,0,0)  =  2n/i '‘r(c/)"  _  ,oxp(  -/let) 


d0' 


II 


M' 
0  ' 


x(  ||  /(/?  >  Vb/ct)  |  /(cos0’()>  Vh/cl)  expl/lfal  +  scc0’(()| 

i  =  l 

x(  |  [  /ftnj  cos0'(([  |  |  sinO'^  1 1’(|  _“|  V-1"--’.  ln>  2)  (53) 
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wherein  it  is  understood  that  the  product  in  braces  is  to  be  omitted  in  the  case  n  =  2. 

Next  we  make  the  integration  variable  transformation 

,v  -  C  cosO'  v 


/  V  —  cost/  . 

v  -*  v  =  arctanl  -  ).  (i  =  0,...,n  —  2) 

1  1  V  C  sin0'  / 


The  nature  of  this  transformation  can  be  best  appreciated  in  terms  of  the  geometry  of  the  triangle 
formed  by  the  three  points  S  ,  S(  + 1  and  S  ,  as  shown  in  Fig.  4:  If  the  line  through  S(I  perpendicular 
to  the  line  through  S(  and  S(  + 1  intersects  the  latter  in  the  point  Tj(  then  v.  is  just  the  angle  between 
S  T  and  S  S  .  , .  It  can  be  seen  from  Fig.  4,  and  it  can  also  be  shown  from  Eq.  (54),  that  as  v  runs 
from  0  to  oo,  the  angle  v  runs  from  -  («/ 2  -  O'  t)  to  rc/2.  Therefore,  this  transformation  renders  the 
integration  domain  bounded;  we  shall  see  shortly  that  it  also  renders  the  integrand  bounded. 

To  calculate  the  Jacobian  of  the  transformation  defined  in  Eq.  (54),  we  begin  by  solving  that 
equation  for  i/(: 


From  this  it  follows  that 


We  also  note  from  Pig.  4  that 


from  which  it  follows  that 


v  =  C  sinfl'  tanv  +  C  cost)'  . 

It  III  l 


—  =  C  sin0’  sec2v  . 
dv  '  •  1 

l 


C  .  .  cosv  =  T  S  =  C  sin0’  , 

l+l  l  IH  l  l1 


CI  +  1  =  C\  sinflVsecv  . 


Now,  a  moment’s  inspection  of  Eqs.  (55)  and  (57)  will  show  that  depends  on  v(  and  C(,  while  Cf  in 
turn  depends  on  C(_  l  and  v(  etc  ,  and  hence  that  o(  depends  on  v|F  v(_ ,, . ..,  v0,  but  not  on  v[  +  j, 

vj  +  2 . vn  This  implies  that  the  Jacobian  determinant  dtuj/dtvj  has  zero  entries  everywhere 

on  one  side  of  the  main  diagonal,  and  therefore  that  the  determinant  is  simply  equal  to  the  product 
of  its  diagonal  elements: 

d(vn,...,v  .,)  n~'t  dv 

-  =  I  —  (58) 

J(v v  )  1  1  dv 

0  ’t-2  i=0  i 

Taken  together,  Eqs.  (56)  -  (58)  imply  that 


"  -  2  dv 


2-^(1  )=  ||  2** 

V  )  '  1  1  '  /  1  1  ()V  I 

O’  ’  >1-2  1=0  1=0  I 
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ft 

V. 

•% 

1 


►  \ 
r.\ 


=  lie  sin"0'  sec'v  , 

iii  i  i 


i  =  0 

n  —  'Z  .  C  £ 

c 


•2  /I -2 


ii -I 

f,2 

1=0 


c 

I 


But  Eq.  (49)  implies  that  C  _ ,  =  w  ,  and  C'0  =  1 ,  so  we  conclude  that 

J/"-2  v  "~2 

— - 2=1  (  f  lsinO’  o-2 ,  =  Ilf. 

a(V-\-2)V,=0  ,=!.  ' 


(59) 


Eq.  (59)  tells  us  that  the  transformation  (o0,...,  ur/  vn  _2)  defined  in  Eq.  (54) 

brings  the  integral  in  Eq.  (53)  into  the  form 


f  'l’n  r  nn  ,  n-  2  ,  ,,n  ,’n  i 

V,o.o,  =  t.#.W W-JM  |o  -r.  j„  |  f  1 1  J,  <»',  |„  I.  *.  | 

0  1  _  1  i 

n  - 1 

x|  [j  /(B  ^  >  V'/Vc/)  j  /(cosfl’()>  V/i/rf)exp|//fc(l  +sec0'())l 


x(  ri  a©,  >)«■»•„(  1 1  )  v — x. 

t=l  i=0 


<n>  2) 


(60) 


wherein  it  is  understood  that  the  product  in  braces  is  to  be  omitted  in  the  case  n  =  2. 

For  n  =  2  the  last  two  factors  in  Eq.  (60)  are  both  unity  (recall  that  C()  =  1 ),  so  the  integrand 
is  clearly  bounded  for  n  =  2.  |The  exponential  involving  sec0'()  in  the  integrand  causes  no 
boundedness  problems,  because  the  last  /-function  in  Eq.  (60)  imposes  the  condition 
bsec0'()  <ct/V<ct.  I  For  3  the  integrand  of  Eq.  (60)  contains  factors  of  C,,...,  f,’(|  any  of  which 
can  be  arbitrarily  large.  However,  the  contribution  of  these  factors  to  the  integrand  is  moderated 
by  the  quantity  V  according  to 


n  -  2 

n  -  2 

IK 

2 

sJ 

II 

ri 

\ 

7 

(o>3) 

i  =  0 

i  =  i 
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The  definition  of  V  in  Eq,  (470  shows  that  V  is  jdst  the  perimeter  of  the  figure  OS,— SO  [see 
Fig.  3|,  and  it  is  obvious  that  this  perimeter  cannot  he  less  than  twice  the  length  of  any  cord  C  ; 
therefore,  CJV s  1/2  for  all  i,  whence 


(  U  ci  )  <  (1/2)”’ 


(n>  3) 


We  conclude  that  the  integrand  in  Eq.  (60)  is  bounded  for  n>3. 

We  collect  below,  in  Eqs.  (62),  the  formulas  through  which  the  various  quantities  in  the 
integrand  of  Eq.  (60)  may  be  calculated  in  terms  of  the  integration  variables.  Eq.  (62a)  follows 
from  Eqs.  (48),  (47a)  and  (51);  Eq.  (62b)  follows  from  Eqs.  (55)  and  (49);  Eq.  (62c)  is  Eq.  (5a);  Eq. 
(62d)  follows  from  Eqs.  (48);  Eq.  (62e)  follows  from  Eq.  (55);  the  components  e(Jt,  c,  y  and  e.tJ  in  Eq. 
(620  are  to  be  calculated  from  C(,  0'  and  tf>' t  according  to  the  formulas  in  Eqs.  ( A8);  Eq.  (62g) 
follows  from  Eqs.  (48);  Eqs.  (62h)  and  (62i)  both  follow  from  the  second  of  Eqs.  (49);  Eqs.  (62j)  and 
(62k)  both  follow  from  Eq.  (47e);  Eq.  (621)  is  Eq.  (470;  and  finally,  Eqs.  (62m)  and  (62n)  follow  from 
the  definition  of  0  as  the  angle  between  ef  ,  and  e(,  together  with  the  fact  that  eri  =  -  a  =  -  C(). 

The  geometric  content  of  the  formulas  below  is  summarized  in  Figs.  3  and  4. 

C0  =  x  sin0’o  +•  z  cosf/()  (62a) 

vn  =  sin0'()  tanv()  +-  cos  O' ^  (62b) 

e()  =  z  (62c) 


C,  =  C.-I  “  Vle<-i 

vt  =  Ct  sin0'  ■  tanv{  +  C(cos0'( 
e  =  xe  +  ye  +  ze 

i  IJC  J  t.y  U 


C  i  =  C  -  v  e 

n - l  n  - 2  n- 1  n -  1 


(n>3:  i  -  1,..  ,  n-  2) 


v  ,  =  C  . 

n  -  1  n  -  I 


>  .  =  c  ,/c  . 

n  -  1  »i  -  l  n  -  1 


_  uoeo 

B  =  B  .  f  c  .  e  . 

I  I  ~  1  1  -  I  /  -  1 


(n  >3;  i  =  1 n  -  1 ) 


n  -  l 

V  -  1  +  \ 


0  =  arccosle  ,-e)  (i  —  1,...,  rt  —  1) 

l  f  -  If  ’  ’ 


(62m) 


0  =  arccos(  -  e  ,  •  C..) 


(62n) 


NWC  TP  6605 


Notice  that  Eqs.  (62d,e,f)  and  Eq.  (62k)  are  not  used  if  n  =  2  The  interrelated,  recursive  structure 
of  the  formulas  for  C.,  u(and  e(  in  Eqs.  (62a  -  i)  would  make  the  derivation  of  explicit  formulas  for 
those  quantities  quite  complicated,  especially  if  n> 3;  fortunately,  explicit  formulas  are  not  needed 
for  computational  methods  that  utilize  a  digital  computer. 

Eqs.  (60),  (62)  and  (A8)  are  quoted  in  Ref.  1  as  Eqs.  (44),  (45)  and  (46),  respectively 

5.  CUBING  THE  INTEGRATION  DOMAIN 

The  integral  expression  in  Eq.  (60)  is  free  of  delta  functions  and  has  a  bounded  integrand 
and  a  bounded  integration  domain.  It  is  therefore  amenable  to  direct  evaluation  bv  standard 
numerical  methods.  However,  many  numerical  methods  arc  easier  to  implement  if  the  integration 
domain  is  the  unit  cube.  In  this  section  we  shall  make  a  final  change  of  integration  variables  to 
bring  Eq.  (60)  into  the  form  of  an  integral  over  the  (3n  -  4)-dimensional  unit  cube 

Each  of  the  integration  variables  fl'(,  v(  and  in  Eq.  (60)  measures  an  angle  that  has  a 
clear  physical  interpretation  in  terms  of  the  geometry  of  the  path  of  an  n-scattered  photon  1  see 
Figs.  3  and  4|.  We  are  now  going  to  change  from  these  integration  variables  to  a  new  set  of 
variables,  p.,  qt  and  wt,  whose  physical  interpretations  are  quite  obscure  but  whose  lower  and 
upper  integration  limits  are  all  0  and  I,  respectively.  The  particular  transformation  that  we  shall 
use  for  this  purpose  also  has  the  convenient  property  that  its  Jacobian  is  a  constant,  this  means 
that,  apart  from  a  different  factor  in  front  of  the  integral,  the  integrand  in  Eq.  (60)  will  be 
unchanged  by  the  transformation  The  actual  derivation  of  our  transformation  equations  uses  a 
special  analytical  technique  in  Monte  Carlo  theory  called  the  "generalized  inversion  generating 
method.”  This  analytical  technique  is  discussed  in  detail  in  Secs.  2-5,  2-6  and  4-6  of  Ref  3  The 
derivation,  although  not  particularly  difficult,  is  moderately  lengthy;  therefore,  we  shall  be 
content  here  simply  to  state  the  result  and  then  verify  that  the  transformation  indeed  has  the 
special  properties  claimed. 

With  k0  defined  by 

Si  s  - v«)/«2,  <63a) 

we  define  the  variables  p{)  and  q{)  so  that 

«'0  =  «H  -  (I  -,c0p(|),/'-'l.  (63b) 

v()  =  «|I/2  -  t/0(  1  “I  <63c) 
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We  also  define  the  variables  pt,  qt  and  wi  for  i  —  1  -  2  so  that 


=  bH  -  p"2), 
v  =  n(l/2  -  qj),1'2),  ‘ 
<t>’  —  2 nw  . 

i 


(n>3;  i=  - 2) 


(63d) 

(63e) 

(63f) 


Consider  Eqs.  (63b-c).  The  first  of  these  two  equations  shows  that,  as  p()  runs  from  0  to  1 ,  W() 
runs  monotonically  from  0  to 

nil  -  < l  - /c{J),/“ I  =  (/»„, 


where  the  equality  follows  upon  substituting  for  *0  its  definition  in  Eq.  (63a).  Eq.  (63e)  shows  that 
for  p0  and  hence  also  0'0  fixed,  as  q0  runs  from  0  to  1,  v0  runs  monotonically  from  n/2  to 

«|  1/2  -  (1  -K0p0)U2\  =  0'n  -  ill 2, 

where  the  equality  follows  upon  substituting  ford  -  /c0p0)I/2  from  Eq.  (63h).  We  conclude  that 
Eqs.  (63b-c)  map  the  unit  square  in  p()q()- space  onto  the  two-dimensional  region  in  0'()v()-space 
defined  bv  the  integration  limits  on  0'()  and  v()  in  Eq  (60;  Since  0'()  is  independent  of  q(),  then  the 
corresponding  Jacobian  of  this  subspace  transformation  is 


rlO'  (IV 

O’  0 

0  0 

^0 

=  |/< <l/2)(l  -Kf}pn)  1/2(;c())||fi(l  -K(Jp0),/2|  =  k()u2/ 2 


Turning  next  to  F^qs.  (63d-fi,  the  first  of  these  equations  shows  that,  as  p;  runs  from  0  to  1. 
O'  runs  monotonically  from  n  to  0.  F^q  (63e)  shows  that  for  pt  and  hence  also  O'  t  fixed,  as  qt  runs 
from  0  to  1 ,  v  runs  monotonically  from  n/2  to 

ri(  1/2  -  p/'2)  =  0\  -  n/2 , 

where  the  equality  follows  upon  substituting  for  p(l/2  from  Flq.  (63d).  And  finally,  Flq  (63f)  shows 
that,  as  u>i  runs  from  0  to  1 ,  runs  monotonically  from  0  to  2 re.  We  conclude  that  Flqs.  ((j.3d-D 
map  the  unit  cube  in  ptqw- space  onto  the  three-dimensional  region  in  O'  \  <p' f- space  defined  by 
the  integration  limits  on  O'  v(  and  ip'  ( in  F’.q  (60)  This  mapping  is  one-to-one  everywhere  except 
on  the  plane  defined  by  p  —  0,  which  is  mapped  onto  the  line  defined  by  O'  t  -  it  and  v(  =  n/2.  Hut 
since  that  plane  and  line  have  zero  volume,  this  lack  of  strict  one-to-oneness  has  no  effect  on  the 
three-dimensional  integrals  of  interest  to  us  here 


For  i  >0  we  note  from  FJ|s  (63d-e)  that  O'  t  is  independent  of  both  t/(  and  u\,  while  v(  is 
independent  of  w .  It  follows  that  the  Jacobian  of  the  transformation  between  0’  x<p' t  space 
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and  p.q-u»-- space  defined  by  Eqs.  (63d-f)  has  zero  entries  everywhere  on  one  side  of  the  main 
diagonal.  Therefore, 

a(0y v*V 

Finally,  since  there  is  no  cross-coupling  between  integration  variables  with  different  index 
values,  the  Jacobian  of  the  full  transformation  is  just  the  product  of  the  subspace  Jacobians. 
Therefore, 


(KJ  av  om 

i  i  i 

d  D  da  dw 


=  [ir  (1/2)  p  ~  1/2|  [n  pl/2|  [2«|  =  a3. 


d^O,VO’0  rVt’^  l’”’®  n-2’Vn-2,<t>  n-2 * 


n  k 


2  t  3v/t  —  2  ^  3/1  —  4  IC,  .  . 

—  (a  )  -  ~  ii  icn.  (64) 

2  2  0 


d(p0,q0,pvqvwl . P„_2.9„_2.«"n_2) 

We  conclude,  then,  that  the  integral  in  Eq.  (60)  transforms  under  Eqs.  (63)  to 


*  (<,0,0)  =  re3*'1  UK  P"c(ct)n  3exp(-/?c()  f  dp  f  dq  \  |  f  dp  [  dq  [  dw 

Jo  J  0  l  ,  _  ,  >  0  J  0  J  0 

.  n-l 

xy  f(B  ^  >  V6/c()  j/(cos6'0>  V,ft/c/)exp|/?/)(l +sec0'{))| 


x(  i  i /*»,»)  cw»-0(  iicjv--* 

1=1  1=0 


(n  ar2)  (65) 


In  this  our  final  expression  for  P  it, 0,0),  it  is  understood  that  the  product  in  braces  is  to  be  omitted 
in  the  case  n  =  2,  and  also  that  the  integrand  is  to  be  evaluated  in  terms  of  the  integration 
variables  through  the  formulas  listed  in  Eqs.  (62)  and  (63)  [see  also  Figs.  3  and  4| 

Eqs.  (63),  (64)  and  (65)  are  quoted  in  Kef.  1  as  Eqs.  (47),  (48)  and  (49),  respectively. 
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FIGURE  I.  Trajectory  of  an  n-Scattered  Photon.  The  photon  leaves  the  origin  O  along  the 
positive  z-axis,  scatters  exactly  n  times  in  the  cloud,  and  then  returns  to  O  at  an  angle  if)  with  the  z- 
axis  The  «th  scattering,  through  polar  angle  0  and  azimuthal  angle  <pt,  occurs  at  point  S(.  The 
vector  from  S(  to  S(  t ,  is  denoted  by  u(  se,u,.  where  ef  is  a  vector  of  unit  length  and  S0  =  Sn  +  1  =  0 
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A9’  * 


7/b,  A 


(a)  n=2 


(b)  n=3 


(c)  n=4 


FIGURE  2.  Geometric  Interpretation  of  the  Relations  Among  the  Principle  Variables  in  Eqs.  (46) 
and  (47)  for  (a)  re  =  2,  (b)  n  =  3,  and  (c)  n  =  4.  The  vector  v(,  with  magnitude  i>(  and  unit  direction  e  , 
represents  the  scaled  path  of  the  photon  between  the  ilh  and  (i  +  l)th  scatterings.  The  angles  0  and 
<pt  are  the  polar  and  azimuthal  angles  of  e;  relative  to  the  polar  direction  e(  _  r  The  vector  OS"  for 
i  —  1,...,  n-  1  is  designated  B  ;  the  vector  ()S^  is  designated  a,  and  has  unit  length  and  polar  angle 
ip.  The  xyz-frame  is  defined  so  that  e0  =  z  with  a  lying  in  the  xz-plane  Notice  that  the  quantity  V 
defined  in  Eq.  (470  is  the  perimeter  of  the  (generally  non-planar)  figure  OSj---S  O. 


■  VvV-*’ 
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vi=ci. 


(a)  n  =  2 


(b)  ii=3 


(c)  «= 4 


FIGURE  3.  Geometric  Interpretation  of  the  Relations  Among  the  Principle  Variables  in  Eqs. 
(53),  (60)  and  (62)  for  (a)  n  =  2,  (b)  n  =  3  and  (c)  n  =  4.  The  angles  6’  and  < p .’  are  the  polar  and 
azimuthal  angles  of  v  =  uei  relative  to  the  polar  direction  C(  =  S(Sa  (t=l,...,n  -  2);  the  other 
quantities  are  as  specified  in  Fig.  2.  (The  vectors  B(  =  OS(  are  still  present,  but  they  are  not  shown 
here  in  order  to  avoid  complicating  the  diagrams.)  Note  that  a  and  tf>  in  Fig.  2  have  here  been 
renamed  Cn  and  6n',  respectively. 
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APPENDIX:  COMPONENTS  OF  THE  VECTORS  e.,  IN  THE  xyz-FRAME 


Let  C  be  a  vector  with  components  Cx,  Cy  and  C2  relative  to  some  coordinate  frame  F  with 
basis  vectors  x,  y  and  z.  Define  frame  F'  as  that  frame  whose  basis  vectors  x’,  y'  and  z'  satisfy 
z'«C  and  y'«z  xC  [see  Fig.  A1 1.  Finally,  let  e  be  a  unit  vector  having  polar  angle  O'  and 
azimuthal  angle  <p’  relative  to  frame  F\  We  want  to  calculate  the  components  ex,  ey  and  e^ofe 
relative  to  frame  F  in  terms  of  the  quantities  Cx,  Cy,  C2,  O'  and  <p' . 


We  begin  by  defining  the  auxiliary  quantities  ct,  cy,  c^  and  c  by 


cxs(y<;-  cy-CrlC.  c2  S  C  /c.  cxv=  (C/  +  C*)»«/C. 


2 +  c  2)112 

Then  the  polar  angle  //  and  azimuthal  angle  4  of  C  in  the  F-frame  are  given  by 


(All 


cos//  =  cz,  sin//  =  cxy  (A2al 

cos*  =  cx  lcfy,  sin*  =  cy/c  .  (A2hl 

Fig.  A1  shows  how  the  angles//  and  *  determine  the  orientation  of  frame  F'  relative  to  frame  F. 
From  the  geometry  of  that  figure,  we  can  see  that  the  projections  of  the  primed  unit  vectors  onto 
the  unprimed  unit  vectors  are  as  follows: 


x’-  x  =  cos// cos*  =  czcx!c  , 
x'-y  =  cos// sin*  =  cgcv/c  , 
x'-z  =  -sin//  =  -  c  ty; 

y'-x  =  -sin*  =  -cy/cxy, 
y’-  y  =  cos*  =  cx  icxy, 
y'-z  =  0; 

z'  -  x  =  sin// cos*  =  cx, 
z'-  y  =  sin//  sin*  =  c  , 
z'-  z  =  cos  (I  =  cy. 


<A3a) 


( A3h) 


( A3c) 


Also,  since  O'  and  <p'  are  defined  as  the  polar  and  azimuthal  angles  of  the  unit  vector  e  relative  to 
frame  F',  then  the  projections  of  e  on  the  primed  unit  vectors  are 

x'>  e  =  sinO'  cos0',  I 

y'-e  =  sinO'  sin</>',  }  (A4) 
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FIGURE  Al.  Relative  Orientation  of  the  xyz-Frame  and  the  *Vz'-Frame.  The  latter  frame  is 
defined  by  z'*C  and  y'^zxC.  The  angles  \i  and  i  are  the  polar  and  azimuthal  angles  of  the 
vector  C  relative  to  the  xyz-framc. 
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Now  we  observe  that  the  x-component  of  e  can  be  calculated  as 

ex  =  x-e  =  Ix'(x'-x)  +  y'(y'-x)  +  z'(z'-x)|*lx'(x'-e)  +  y'(y’-e)-t-z’(z'-e)|, 


or 


(A5a) 


ex  =  (x'-x)(x'*e)  +  (y’-x)(y'-e)  +  (x’>x)(z’-e). 

Similarly,  the  y-  and  z-components  of  e  can  be  calculated  as 

(A5b) 

ex  -  (x'*z)(x'*e)  +  (y'-z)(y'-e)  +  (z'*z)(z'-e).  (A5c) 

If  we  now  substitute  Gqs.  (A3)  and  ( A4)  into  Eqs.  ( A5),  we  find  that  the  resulting  equations  for  e  , 


ey  =  (x'-y)(x'-e)  +  (y'*y)(y'-e)  +  (z'-y)(z'-e), 


e  and  e  can  be  written  in  matrix  form  as 
y  * 


c  c  Ic  -c  /e 

z  x  xy  y  xy 


c*cy'c*y 


Cx,Cxy 


-C 


xy 


y 

c. 


sin0'  cos <p' 
sin0’  sin#' 
cos  O' 


(A6) 


Together,  Eqs.  (A6)  and  (Al)  give  the  F-frame  components  of  e  in  terms  of  the  F-frame 
components  of  C  and  the  angles  0’  and  <t>'.  But  notice  that,  if  C  should  happen  to  coincide  with  z, 
then  c  =c  =c  =  0,  and  the  four  quotient  elements  in  the  above  3X3  matrix  become 

x  y  xy 

indeterminate.  Since  in  that  special  case  we  will  have  z'  =  z,  then  we  may  as  well  take  x’  =  x  and 


y'  =  y.  Therefore,  if  cxy  =  0  then  we  shall  simply  take  the  3x3  matrix  in  Eq.  ( A6)  to  be  the  unit 


matrix  (with  l’s  along  the  main  diagonal  and  0’s  elsewhere). 

The  foregoing  result  is  actually  used  in  two  different  ways  in  the  text  to  calculate  the  xyz- 
frame  components  of  the  unit  vectors  e-.  The  first  way  is  in  connection  with  Eqs.  (5)  and  (6). 
There,  e(  (i  =  1 ,...,«)  is  stipulated  to  have  polar  angle  0(  and  azimuthal  angle  <Pi  in  the  coordinate 
frame  whose  z-axis  points  along  e(  (  and  whose  y-axis  points  along  z  x  e(  _  , .  Therefore,  if  eiJC,  e(  y 
and  are  the  components  of  e(  in  the  xyz-frame,  and  if  eiJcy  =  (e ,^2  +  e ,  >2)l/2,  then  above  result 
implies  the  recursion  relation 


ij t 


e  ,  e  ,  le  ,  -  e  ,  le  , 

i  -  1  ^  i  -  1  ,.f  i  -  I  ,x>  i  -  l  ,y  <  -  I  jcy 


I-1.J 


e  .  e  .  !e  . 

t  -  1  j.  i  -  1 ,  v  i-l  jc  y 


VuVtjcy 


c‘-i.y 


•sj 

where,  from  Eq.  (5a), 


i-l  jcy 


i-U 


sin0(  cos  <pt 

sin0(  sin#( 

cos  0( 

, (t=l,..., it)  (A7a) 


Ojc 


=  °-  'o.v^0’  \x=1’ 


( A7b) 
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and  where  it  is  understood  that,  if  eiljty=0,  then  the  3X3  matrix  in  Eq.  (A7a)  is  to  be  taken  to  be 
the  unit  matrix. 


Eqs.  (A7)  show  how  the  xy2-frame  components  of  e.  may  be  calculated  from  the  xyz-frame 
components  of  e •  _ ,  and  the  angular  variables  8i  and  #,•  These  relations  are  required  for  a 
complete  interpretation  of  the  "early”  formulas  in  our  analysis  [specifically,  Eqs.  (1)  through  (46)  |. 
However,  our  final  formulas  for  P^U.0,0)  for  n  s  3  [Eqs.  (60)  and  (65)1  have  as  their  integration 
variables,  not  the  angles  0(  and  but  the  angles  0.'  and  <pl  (t  =  l,...,n-2).  These  primed  angles 
are  defined  [cf.  Eqs.  (52)1  as  the  polar  and  azimuthal  angles  of  e.  in  the  frame  whose  z-axis  points 
along  C(  and  whose  y-axis  points  along  z  x  C  ■,  where  C(  is  the  vector  defined  in  Eq.  (48).  It  follows 
from  the  foregoing  analysis  that,  in  Eqs.  (62),  the  xy2-frame  components  of  the  vectors  elt...,  e„  _2 
are  to  be  calculated  according  to  the  following  formula: 


where 


~eiJC 

c  c  lc  -c  lc 

u ijc  ijcy  i.y  Ucy 

Cijc 

sinQ^  cos#,' 

et.y 

= 

c  c  lc  c  lc 

tj  I.y  Ijcy  tjc  Ijcy 

c‘.y 

sin6('  sin#/ 

-c  0 

uy 

ci;_ 

cos  O' 

i 

,  U=  1 . n-2)  (A8a) 


c  ■C.  1C,  c  1C,  c  *  C  1C,  c-  ■  (C  Z  +  C  ‘)irz/C, 

jjt  ijt  i’  i.y  i’  m  i *  i’  ucy  ijc  i.y  i’ 


2-*-/’  2'1'2"'’  (A8b) 


and  where  it  is  understood  that  if  c,  „„  =  0  then  the  3x3  matrix  in  Eq.  ( A8a)  is  to  be  taken  to  be  the 
unit  matrix. 


Eqs.  (A8)  are  quoted  in  Ref  1  as  Eqs.  (46). 
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